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Victor smoothing splines on the sphere are deftned, 

Theoretical properties are briefly alluded to. An approach to choosing 
the appropriate Hilbert space norms to use In a specific meteorological 
application Is described and Justified via a duality theorem. Numerical 
procedures for computing the splines as well as the cross validation 
estimate of two smoothing parameters are given. A Monte Carlo study 
Is described which suggests the accuracy with which upper air vorticity 
and divergence can be estimated using measured wind vectors from the 
North American radlosoncte network. 


To appear In **Mult1 variate Approximation Theory", Vol. 2, proceedings 
of a conference at Oberwolfach February 7-13, 19^. H. Schempp and 
K. Zeller, eds. Birkhauser Verlag, Basel -Boston-Stuttgart. 



1. Introduction 


A thoory of splint functions on tht sphtrt is rtpidly btinp drrtloptd. 
stt UAHBA (1981i), FREEDEN (1981i.b). SHURE. PARKER AND BACKUS (1981). 

Dr. FREEDEN will bt rtportin^ on soM of his rtsults tistwhtrt in this 
volunt. Much of tht rich thtory surrounding univtrittt spUnts tnd 
thin platt splints cittrly is txttndthit to tht thtory of splints on tht 
sphtrt* via tht ust of reproducing ktmtis* n-widths* ttc. In particular 
convtrgtnct ratts for smoothing splints on tht sphtrt can bt obtaintd 
fran tht known ratt of dtcay of tht tigtnvaluts of tht rtitvant reproducing 
kernels, sat a.g. NICCNELUand NAHBA (1981), UAHBA (1977), UTRERAS (1981). 

In this paper wt propose a notion of vector splints on the sphtrt. It 
is clear that interesting approx inat ion theoretical properties of these 
splines can be obtained. However, in this paper our focus will bt on the 
solution of certain practical problems which must bt solved so that these 
splines may be usefully applied to the analysis of neteorological data 
from the upper air radiosonde network. 

For the purpose of numerical weather prediction the global radiosonde 
(weather balloon) network takes measurements every 12 hours of the 
horizontal wind velocity vectors and other variables, at 9 standardized 
vertical levels. From this data it is desired to estimate the horizontal 
wind field and its vorticity and divergence (and other variables) at 
a regular grid of points, for each level. These estimates on a grid art 
then merged with estimates of the same variables cn the same grid, which 
have been obtained from a forecast, to provide an estimate of the present 
state of the atmosphere. This state estimate is then used as the Initial 
conditions co a manerical integration scheme which Integrates a set of 
differential equations describing the dynamics of the atmosphere, to 
provide a new forecast. Numerical weather forecasts can be quite sensitive 
to errors in the vorticity and divergence in the initial wind fields. 
CNifortunately, horizontal wind vectors at, for example the 500 millibar 
height, of the order of a few tens of meters per second, are measured 
with an error standard deviation in each component of the order of 
2*4 meters per second. Thus, it is not a trivial matter to obtain accurate 
information concerning the vorticity and divergence from this data, even 
In areas such as the continental U.S. where the radiosonde network Is 
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rtUtlvtly dtntt, Nt telitvt tMt ^ «pproprtitt of Hii voder 

SNOthlng sellnes d propose hove tlio potontlol for ddiip tMs roloMvoly 
well. 

Spooklnp Intuitively, tlio vector smo^Ii^ tpllnet ee prepeee frill 
behive like 1«r peso filters. Xn tin tpllnet nr propose tliere trill be . 
two regulsHtstlon or seoothtng ptroMters tj cVtsee snd two (eets of) 
"shspe" perdueters. The first seootlilnp permeter to bo cbeten. esy 
be thought of ss governing the overell hslf power poirt of the low 
piss filter. The second pareneter governs the relative distribution of 
power between vortictty and divergence In the estlMte. The choice of 
the two sets of shape parmeters correspond to the choice of Hilbert 
space norms, but In an Important practical sense they govern the rates 
of decay of the energy sp^trum of the solution, one "shape" for verticity. 
and one "shape" for divergence. It Is well known from the theory and 
practice of 111 posed problems that the appropriate choice of certain of 
these parameters can affect the practical usefulness of the result. 

In this paper we propose the use of generalized cross-validation (GCV) 
for choosing the two smoothing parameters. GCV can also be used to 
choose a snail number of "shape" parameters (see CRAVEN and MANBA (1979). 
WAHBA and UENDELBER6ER (1980)). However. In this paper we show how 
historical meteorological data can be used to choose the "shape" 
parameters, or Hilbert space norm. We discuss some numerical methods, 
and we describe the results of some numerical experiments on synthetic 
data which mimics actual 500 millibar horizontal wind fields over the U.S. 

In our experiments we have observed that the accuracy In estimating both 
vorticity and divergence can be quite sensitive to the relative distribution 
of power allocated between them, (choice of second smoothing parameter) but 
that GCV can be quite effective In estimating the correct relative power 
distribution. 

For the meteorological experts In the audience we remaHc that 
estimating the present state of the atmosphere from current data Is not 
exactly the same problem as ostlmatlng the state of the atmospiMre from a 
combination of present data and a forecast of the present. This Is so 
because a data only estimate needs to take account of properties of the 
atmosphere and measurement system while a data plus forecast estimate needs 
to take Into account the relative error of the data and the forecast. In 
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thit p«ptr «• irt stu«(ytn 9 thi 4«tt cnly preh1«i, Hewtv«r, m btlltM 
this class of ttchnlquos can bt axtandad to tha data plus fortes probl«i 
and hope to do that In a subsaquant papar. 

In Saetton 2 wa dafina tha vactor saoothlng spllnas. In Sactioa 3 
wa discuss tha choica of Hllbart spaca nonas. In Sactlon 4 wa dascriba 
nuRarlcal nathods and tha cross validation astliaata of tha'saoothlng 
paraiRetar end In Saetton 5 wa dascriba a Nonta Carlo test of tha Mthod. 


2. Helaholz Thaoraa and Tha Definition of Vactor Snoothl no Splines 

Wa let P be a point on the sphere S. P ■ (X,4), where X ■ longitude 
(0<X<2w) and ^ ■ latitude (-2<4<|). V ■ (U,Y) Is a (sufficiently racier) 
horizontal vector field on the sphere, wl«re U(P) is the eastward component 
and V(P) Is the northward component at P. 

The vorticity c and the divergence D of V are given by 


where a Is the radius of the sphere. Then there exists (by HelmoUz 
Theorem) two functions f(P) and o(P), PeS, called the stream function and 
the velocity potential respectively, with the following properties: 


U 


1{. H + 

a^ cos(^ ax^ 


V 


1/ 1 ay Wv 
a'cos(j> ax 


(2.3a) 


D ■ A4 


(2.3b) 


where a Is the (horizontal) Laplaclan on the sphere 


1 r 1 
Af ■ - 


1 


^XX * cos4> 


<cos*f^)^]. 


a^ cos*4i 

f and are uniquely determined up to a constant (which we will take to be 
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determined by /'F(P)dP - /«(P)dP -0, Ue are interested In defining Hilbert 
S S 

spaces of vector fields whose divergence and vorticity exists pointwise. Me 
will do this as follows. Let f(P} be a square Integrable ^ncti(m on the 
sphere which integrates to 0. Then f has an expansion in the normalized 
spherical harmwiics 


00 

I I Vi' 

1-1 s-t “ ’■ 


where 


9j^cossXP^s1n<^) 0<s<il 

» t » 1,2,... 

9^sinsXPj|^'(s1n(})) -t<s<0, 



and the Fourier Bessel coefficients fj^ are given by 
fjj . /f(P)Yj*(P)dP 

with 

/f"(P)dP . I f‘. 

Ifi 

s 

Now Y are the eigenfunctions of the Laplacian 

X 

AYj^^ = -t(t+l)Y^^ . 

Thus 

Af = 

£ • s 

Let i » 1,2,..., s * be a set of nonnegative numbers with 

\ » max X and 

^ S=-t,...2- ^ 
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(2.5) 

I *■ 

Using the addition formula for spherical harmonics 



where y Is the angle between P and Q, the Cauchy-Schwartz Inequality and 
the fact that P|^(7) • 1 gives 


|Af(P)|-l I ^1 

X-.5 


I Sm.l 


i(?irZ‘Mt+i)M2ui)i )’/^( Z 

I ^ t,s ^As 

f«5* 

Thus {X.J satisfying (2.5) and I -r — <« Imply that Af(P) Is well defined 

£,S ^£S 

and finite for all P. 

Let H be the collection of all pairs (f,*) on the sphere which Integrate 
to zero, are squar*e Integrable and 

• Z Z 

£■1 ^ ^ 

- I I • WP)Yj^(P)dP 

£»1 S«-£''as' ' *' 

where and {Xj^^(2)} are sequences satisfying 

ao 

y tM£+l)'(2£+l)maxA, (i)«». i = 1,2. 

£“1 S ^ 


H is clearly a Hilbert space with square norm 
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for any fixed d > 0 and both fflembers of each pair possess Laplaclans every- 
where. It Is easy to show that If then 


J(f) 


/(A"‘^f)*dP 

/{ * JdP 

sin"* ♦ 


m even 
m odd. 


( 2 . 6 ) 


If Xj^g ■ Ct(t+1)3*", them m > 3 guarantees the pointwise existence of the 
Laplaclan. 

The observations are assumed to be of the form 

- U(P^) + - V(P^) + e][. 1 • 1.2.. ...n (2.7) 


where (U(P^) .V(P^)) is the true (wind) vector at P^ and e^.e? are 
measurement errors. We propose estimating the stream functl£) and velocity 
potential (f .*) associated with U and V by finding (T.I^)eH to minimize 






1 n ^ 

* a34i''r"'1 

+ xCo^Ch*) + -J2(*)] 


(2.8) 


Note that In the residual sum of squares above, U(P^) and V(P^) are 
expressed In terms of y and «I> via (2.3a), A unique ralnimizer ('{'^^ ^) 

exists for each A > 0, 6 > 0 and the resulting wind field (U, -,V. .) 

A ^0 A 1 0 

constructed from (f, ,) may be termed a vector spline field. Its 

A ft A •O 

vorticity and divergence will be given by 'a.6 ■ 

(Obviously, Interpolating splines can be defined as mlnimizers of 

Jl('l') + subject to the Interpolating conditions, we will not discuss 

these further.) Using WAHBA (1981a) or FREEOEN (1981a) It Is straightforward 
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to write an explicit (Infinite serlesl expression for (U^^ 


3. On the Choice of and 


Let X 


ts 


I I a.[(-i)(£+l)]^r2 and suppose that X,. > 0 for 
j.Q J 


l ■ 1,2,..., s - -£ £. It Is not hard to see that 


0(f) 


00 S, 

'J, ^ 

A«1 


lA 


/( " 
j-0 ^ 


(3.1) 


SO that the choice of the X^^ can then be reduced to the choice of m and 
the (oj)* (If ” 0 for one or more s, the minimization problem can 
be handled by the methods described in KIMELDORF and WAHSA (1971), see 
also FREEDEN (1981a). In principle m and possibly (with 
can be chosen by cross validation (see WAHBA and UENDELBERGER (1980) but 
It is undesireable to attempt to choose too many of these parameters from 
the data, see WAHBA (1981c). 

in this section we will use the duality theorem which relates 
smoothing by splines to Bayesian est 1m at Ion /Weiner filtering on stochastic 
processes to suggest how the J's may be chosen based on historical 
meteorological data. 

To give the duality theorem we need some background, which we will give 
In a univariate context. 

Let X(P). PeS be a (univariate) zero mean Gaussian stochastic process 
on the sphere with covariance R(P,Q) defined by 

R(P,Q) « EX(P)X(Q), 

where E Is mathematical expectation. Following PARZEN (1961), CRAMER and 
LEADBETTER (1967) we can define the Hilbert space X spanned by X(P), PeS, 
as all finite linear combinations of random variables (r.v.'s) of the form 

and their quadratic mean (q.m.) limits. (A sequence Z.j,Z2 ,..., of r.v.'s 
has a q.m. Mmlt If 11m E(Z.-Z^)^ " 0). The Inner product In X Is 

z ,m-K» 
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<X(P).X(Q)> - EXCP)X(Q) ■ R(P,Q1, «nd ts Mtinded by Itnwity to all r.v.'s 
"k 

of the form Z|^ - J thotr q.m. limits. For example . letting 

L be a linear functional, the r.v, LX ■ ^(Pq) *^111 be In X If the se<Hience 
of r.v. 's 

X(XQ,(^Q+h|^)“X{XQ »^q) 

h KT 


(3.3) 


has a q.m. limit, as hi^-^), where (Xq,^q) ■ Pq. Then, It Is not hard to 
show that the sequence {Zj^} will have a q.m. limit Z ■ "^(Pq) and only If 

Is well defined and finite. Then the quantity In (3.4) Is equal to 
E(-^X(Pq))^» and furthermore 

E(^X(P„))X(Q) .^R(P(,,Q) . 

More generally, let be the reproducing kernel Hilbert space with 

reproducing kernel R. Then each random variable of the form Z ■ LX can be 
Identified with the bounded linear functional L on and vice versa. The 
argument Is as follows. If Z ■ LX Is a r.v. In X It can be shown that 
EZX(Q) ■ L^pjR(P,Q) ■ n(Q)t say, where L^pj means the linear functional L 

applied to R considered as a function of P. However, by the properties 
of reproducing kernels, It can be shown that n(*) 1$ the representer of L 
In R, that Is Lf ■ <n,f>p, where Is the Inner product In Hp. We 

are now ready to state the 

Duality Theorem (KIMELOORF and WAHBA (1970)). 

Let X(P), P&S be a zero mean Gaussian stochastic process with 
covariance bR(P,Q), and let Hp be the reproducing kernel Hilbert space 
with reproducing kernel Kp. Let 

■ L^X + , 1 ■ 1 ,2... . ,n. 
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where 1 • 1,2,... ere n r.v.'s In X, end the c^,...,c^ art 
Independent, 0 nean Gaussian r.v.'s. Independent of X(P), PeS, with cqmmxi 
variance o^. Then the conditional expectat1<wi of X(Q), given 
, 1 ■ l,2,...,n, 

E{X(Q)|Y^ ■ y^. 1 • 1.2,... ,n} (3.5) 

Is given by f^(Q)» where fp^(>) Is the solution to the nininizatlon problea: 
Find feHj^ to minimize 

and X ■ a*/nb. 

Proof: See KIMELDORF and WAHBA (1970,1971), WAH6A (1978). However, the 
proof proceeds by direct calculation of fj^(Q) and by using the facts that 
E(L^X)X(Q) » n^(Q), where 

Now let f be seme atmospheric variable of Interest. Me will proceed 
as though the different realizations of f were sample functions from a 
zero mean Gaussian stochastic process with covariance R(*,0« If repeated 
(Independent!) observations on f were available, then various properties 
of R could be estimated from this data. He will discuss both "frequency 
domain" and "space domain" methods for doing this. Using ^he properties 
of reproducing kernel spaces (see, e.g. HASHED and HAHBA (1974)) It is 

f * 

not hard to show that If J(f) • J Is the norm on a reproducing kernel 

t,s '^is 

space H, then the reproducing kernel R for H is given by 

R(P.Q) - I t3.6) 

t ,s 

To simplify the discussion, in this paper we are considering only R's 
whose eigenfunctions are the spherical harmonics. (Other eigenfunctions, 
i.e. those associated with Laplace's tidal equations, may well be reasonable 
in certain meteorological applications, sec WAHBA (1981b)). 
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If a stochastic proctss X(P), Pc$* has covarlanct 
t,s 

then X nay be nodal ltd as a r«tdon linear conblnatlon of the spherical 
hrmonlcs CKarhumen-Loeve expansion 1 

where the X^^ are randan variables with 

“ ^is* • A's', ■ 0, Is f t‘s'. 

(To see this, compute EX(P)XCQ) from (3.7) to obtain (3.6).) We have 

Xj^ • /X(P)Yj^^P) 


and 


\i ■ E*j| * E(/X(P)Yj*(P)dP)^ 


If K Independent observetlons, of a meteorologlcel variable of 

interest are available, this suggests choosing based on estimates 




2 


where the sample Fourler-Bessel coefficients f^^, k • 1,2,...,K are given 
by 


f^ • /f‘^(P)Yj^^(P)dP, 


Figure 1 gives a plot of l*ebruary 1974 monthly 
averages of some atmospheric mean square sample Fourier Bessel Coefficients 
collected by STANF'^RD (1979) fnm Channels 2 and 4 of the Radiometer 
on NIMBUS-5. The radiation received by Channels 2 and 4 respectively can be 
used (crudely) to infer the temperature T(P), P&S In the upper and lower 
stratosphere, respectively below the satellite, By piecing together data 
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k Ic & 

from several orbits, (approximations to) ■ /T (P)Yjj^ (P)dP can be 
obtained. STANFORD has computed monthly mean square values » 



What has actually been plotted In Figure 1 Is the "TEMP SPECTRAL POWER" 
defined as 


TEMP SPECTRAL POWER U) 


1 5 . 
^j-1 




The energy spectrum In Temperature fields Is related to the energy spectrum 
of other meteorological variables, I.e. wind and geopotentlal . We are not 
concerned here with the exact details of these pictures bu. rather that 
sequences can be fitted to this kind of data to provide meteorologically 

reasonable Hilbert Space norms. See KASAHARA (1976) for some plots of 
sample Fourier-Bessel coefficients with respect to the eigenfunctions of 
Laplace's Tidal equations for wind and geopotential. Figure 2 gives a 

plot of an idealized sequence t « 1,2 where was obtained 

by fitting (by an ad hoc procedure), a function of the form 

• 1 1 

to some of the data behind Figure 1. If Xj^^ does not depend on s, 

^ts * covariance 


R(P.Q) “ I ®(P)Y ®(q) 

t,S 

reduces by the addition formula for spherical harmonics, to 

R(P.Q) ' X? I (2U1)X, P (cosy(P,Q)), 

where y(P»Q) the angle between P and Q. Figure 3 gives the function 
p(y) defined by 

CD cet 

p(y) * ^ (ZUDXjPjCcosy)/ I (2ul)x P.(cosQ) 

Jl=l ^ i*l ^ ^ 
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whlch Is associated with the of Figure 2, Figure 4 gives an estiMte 
for p(y) for f(P) « the 500 millibar (geopotenttal) height obtained by 
JULIAN and THIE3AJX (1975) fron sample covariances from data from a network 
of 51 North American weather stations for the winters of 1966 and 1967. In 
estimating p(y)» *n Isotropic covariance function was assumed. The 
purpose of providing Figures 1 and 4 here Is to convince the reader that 
historical collected w* collectable meteorological data may be used to 
choose the norm on K* although the particular data sets exhibited here 
may or may not be the most appropriate. In the numerical experiments 
to be described we have taken the {Xj^j(l)}and (Xj^^(2)} both as In Figure 2. 


4. Numer'^cal Methods. The Generalized Cross»Val1dat1on Estimates of 
X and S . 

Given X,6, ((U^.V^)}. an approximate 

mlnimizer (f,*) of (2.8) can be obtained In the form 




4 > 


N I 




Jl»l s«-t 




a-1 s»-a 


(4.1) 

(4.2) 


where N Is sufficiently large. For other numerical approaches to the 
minimization of (2.8) see WAHBA (1980,1981a), WENDELBERGER (1982). Let 
N I 

N » I 1 “ N^-1 and renumber the Indices (A,s), s ■ -A,..., A, 

A«1 s»-A 


A « as 1,2,...,N. Let X. be the nxN matrix with (1,As)th entry 

<♦> 




and be the nxN matrix with (1,As)th entry 

1 — ] i-Y ^ ( p ) 

a cos<()^ Z\ A ' V 


and let X be the 2nx2N matrix 
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0/f 




X 




Let D^ be the 2Nx2N matrix 


(4.3) 


D 





(4.4) 


where Is the NxN diagonal matrix with ts»tSth entry 1 ■ 1t2. 
Letting z ■ (U^ •«.• »U^*V^ (...iV^)* y ■ (o^ i. .* ». .. i0^) » It Is seen 
by substituting (4.1) Into (2.8) that we have to find y which minimizes 


1 

n 


2 - xvir + ^v'D*\. 


The mlnimizer Is 


Y ■ (X'X+nxO"')‘’*'i 


(4.5) 


By the use of (2.3a) and (4.3), It follows that the estimated wind field 
x.Vy x) at the data points satisfies 

AfO AfO 


’x.aCi' 


A(X)2 


(4.6) 


''x.aC’n’ 

where A(x) Is the 2nx2n "Influence" matrix 
A(X) « X(X'X+nXD"b'^X' 


The generalized cross validation (GCV) estimate of (x,6) Is the mlnimizer 
of the cross validation function V(X,fi) defined by 


1||(i-a(x,6))2|p 
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Thts method for esttmattng smoothtng parameters tn regularization problems 
was proposed In CRAVEN and WAKBA (1979), GOLUB. HEATH and um/< (1979) 
and UAHBA (1977b) and Its numerical and theoretical properties have been 
studied In various places, see for example UTRERAS (1981), We only note here 
the useful property of the 6CV estimate of X and 6. Let the predictive 
mean square error R(X>6). when X and 6 are used be defined by 

• n- 

where U(P.). V(P,) is the true (but unknown) wind vector, and suppose the 
' ^ U y 

measurement errors and are Independent Identically distributed 
zero mean normally distributed random variables. Then under rather 

A ^ 

general conditions, for large n the minimi zer (X,6) of V(X,(5) provides 

a good estimate of the minimi zer of R(X,5). V Is not guaranteed to 

have a unique, or even a finite mlnimizer. Practical difficulties In 

minimizing V though posrible appear to be moderately rare when the 

assumptions are reasonably well satisfied. Various diagnostic tools are 

available In troublesome cases and will be discussed elsewhere. 

The numerical experiment reported In Section 4 was performed on the 

Amdahl at Goddard Space Flight Center, with 2n • 228, N - 15, 2N ■ 448. 

1 /2 

We outline the calculations used. Let W^ ■ XD , and let the singular 
value decomposition (SVD) of W^ be 

W^ • UDj^r (4.S) 


where UU' ■ U*U » ^2nx2n * * diagonal matrix with entries 

bg^. U,{b^} and V are computed using LINPACK. Letting 



then 
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V(X.6) 


, 2n 1 2n 


b.* 2 




, 2n 

2n I 


b. ' 2 


v=l 




b^*+nX 


\ 


“2n 

bjHnX 


M. 


(4.10) 


(4.11) 


For fixed 6,X(6), the nlnlmlzer of (4.10), Is easily found by a global 

A 

search in Increments of logX. Then V(X(5),6) was plotted for 8 values of 
5 chosen In powers of 1/6, and the minimum was readily evident. No doubt 
more efficient and automatic search procedures can be found. 

For large n, N, and poorly conditioned, computing the SVD 
can be expensive, or It can fall to converge In a reasonable time. Some 
shortcut methods which alleviate this problem somewhat and use less storage 
have been developed. (BATES and tfAHBA, (1382) In preparation.) 


5. A realistic Monte Carlo test of the y.ethod 


A number of techniques for estimating divergence of the upper 
atmosphere from radiosonde data have been proposed In the atmospheric 
sciences literature. For example, see SCWIIDT and JOHNSON (1972). In an 
attempt to determine how well the proposed method might work In practice 
a Monte Carlo experiment simulating realistic measured wind data from 
"model" stream functions and velocity potentials has been coded, and 
various experiments run. We describe one such experiment. 

We obtained a model streamfunction and velocity potential of the form 
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by choosing and as nonnally distributed pseudo-randon nunbers with 
mean 0 and variances * ^ts 9^^®” Figure 2. and Cg 

were scale factors chosen so that the simulated c ■ AV and 0 ■ A« had 
magnitudes typical of real atmospheres. 

(/C*dP)^^ - 6 X 10“^/sec.. (/0*dP)^^^ - 1 x 10'® /sec. Hodel wind vectors 
(U(P^) ,V(P^)) were computed from the model (f,*) of (5.1) for {P^} 
corresponding to n ■ 114 North American radiosonde statiwis. The data 

z - (Up....U^, V^). where - U(P^) + - V(P^) + e/. 

were constructed by adding the measurement errors as normally 

distributed pseudo random numbers with mean 0 and standard deviation 
a » 2.5 meters/sec., a realistic value for the measurement error 
standard deviation. Since the ability to estimate divergence will 
depend on the signal to noise ratio. It Is necessary that the values 
of "signal" and "noise" be chosen realistically. The results reported 
here can be expected to be rosier than that obtainable In practice, 
however, primarily to the extent that wave numbers £ > N occur In practice 
but are not simulated here, and (secondarily) because In practice and 
cannot be so precisely matched to the "truth" as they are In this 
experiment. 

Figure 5 shows the simulated wind vectors. Figure 6 shows the 
estimate of the true wind field, plotted on a 5® x 5 ® grid In latitude 
and longitude. Figures 7 and 8 show the model and estimated vorticity 

^ A 

and divergence, respectively. Figure 9 shows V(X(6),6) and R(X(6),6), 

A 

(of 4.8) plotted as a function of 6. In Figures 6-86- 1/36 was used. 

A 

It can be seen that the mlnimizer of V(X(6),5) was a good estimate of 
the mlnimizer of R(X(6),6). Figure 10 gives M$E(c-^ ) and NSE(0ys ) 

and their sun, where X(6),6 X(6),6 




The {P|^} constitute a regular grid inside the U.S. It can be seen from 

Figure 10 that if 6 is taken as too small Cite, divergence Is suppressed), 
then the mean square error in the estimated vorticity becomes large. 




Figure 7: Model and Estimated Figure 8: Model and Estimated 

Vorticity, X 10“Vsec. Divergence, X 10“®/sec 
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An estimate a* for the variance of the measurement error Is fivotlable n 


a»(X.6) 


Tr(I-A(X,6)) 


since the numerator Is the residual sum of squares and the denominator 

A 

Is the equivalent degrees of freedom for error. In this example a was 
2.58m/sec.« very close to the "true" value of 2.5 meters/sec. In those 
occasional sticky cases encountered In practice where V(X,6) has multiple 
minima. If the order of magnitude of a Is known apriorl, the exmalnatlon 
of a can usually be used tc resolve ambiguity. See NAI^A (1981d), 
WENDELBER6ER (1982). Bayesian confidence Intervals are also available 
for these estimates, see Wahba (1981d). 

We have concluded that this approach has much promise for applications. 
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